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Abstract 

The e A method for predicting transition onset is an am- 
plitude ratio criterion that is on the verge of full matu- 
ration for three-dimensional, compressible, real gas flows. 
Many of the components for a more sophisticated, absolute 
amplitude criterion are now emerging: receptivity theory, 
secondary instability theory, parabolized stability equa- 
tions approaches, direct numerical simulation and large- 
eddy simulation. This paper will provide a description of 
each of these new theoretical tools and provide indications 
of their current status. 

1 Introduction 


Robust tools for predicting the location of the onset of 
transition in boundary layers on aerospace vehicles have 
obvious technological importance. For practical engineer- 
ing purposes one desires a prediction tool which is quan- 
tatively and not just qualitatively correct: the issue is not 
whether transition occurs but where. At present transition 
prediction tends to be based on simple correlations, such 
as crossflow Reynolds number or Re$/M\ modified one- 
or two-equation turbulence models which seek to trans- 
late the freestream turbulence level into computations of 
laminar- transitional-turbulent flow; and linear stability 
theory. 

The pioneering work of Smith k Gamberoni [1] and Van 
Ingen [2] introduced an empirical method for estimating 
the location of transition onset based on an amplitude ra- 
tio criterion. This tool has come to be known as the 
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method or the AT-factor method. The next level of tran- 
sition prediction methodology is likely to involve an ab- 
solute amplitude criterion. In order to achieve this goal, 
many more physical effects must be taken into account 
and more analysis tools must be utilized. 

The seeds of transition are the disturbance environment 
in which the vehicle operates. Transition is born by the 
receptivity process in which the background disturbances 
are incorporated within the boundary layer as linear insta- 
bility waves. It is nurtured by a relatively long region of 
linear instability growth. Once secondary instability ef- 
fects develop, they induce rapid growth, and at a suffi- 
ciently high amplitude the nonlinear regimes are entered 
and transition commences, as signified by, say, the skin 
friction or heat transfer rise. 

Characterization of the disturbance background is a pre- 
requisite for an absolute amplitude criterion. Both ampli- 
tude and spectral information are required. This is nec- 
essarily an experimental task and will not be addressed 
further in this paper. 

The theoretical and computational tools which can now 
be brought to bear on the transition problem include 
Receptivity Theory, Linear Stability Theory (LST), Sec- 
ondary Instability Theory (SIT), Parabolized Stability 
Equations approaches (PSE), Direct Numerical Simula- 
tion (DNS), and Large-Eddy Simulation (LES). 

Receptivity theory is a very active area of current re- 
search. The essential problem is that the freestream dis- 
turbances often have much longer length scales than the 
instability waves in boundary layers. Therefore, the incor- 
poration of background disturbances into boundary-layer 
instability waves requires a wavelength conversion mecha- 
nism. A variety of linear and asymptotic techniques have 
been applied to this problem. We refer the reader to [3]- [8] 
for some recent work in this field. The latter two articles 
are particularly concerned with compressible flow. 

This paper furnishes a brief description of the remainder 


1 


of these tools and provides examples of recent work. We 
shall limit ourselves to illustrating these methods for su- 
personic flows, and shall highlight some recent results from 
the theoretical and computational transition program at 
the NASA Langley Research Center. 


Although results will be presented here for flat plates, 
cylinders and cones, the equations and notation will be 
given just for the flat plate. See [9] for the appropriate 
equations for the more general situation. 


2 Governing Equations 


3 Linear Stability Theory 


The starting point for these analysis t ools is, of course, the 
compressible Navier-Stokes equations. In dimensionless 
form the equations for a thermally and calorically perfect 
gas are 
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is the viscous dissipation. The Reynolds number is de- 
noted by Re , the Prandtl number by Pr, the Mach num- 
ber at the boundary-layer edge by A/ e , and the ratio 
of specific heats by 7 . (For all the examples in this 
paper 7 = 1.4.) In these equations p is the density, 
u = (tii , u 2 ,u 3 ) = (ti, v y w) the velocity, p the pressure, 
and T the temperature. We shall denote the solution vec- 
tor by q = (p, u, t>, w,p). The coordinate system is chosen 
so that x = (xi,x 2 ,x 3 ) = where x, y, and 2 are 

the streamwise, spanwise, and wall-normal coordinate di- 
rections, respectively. 

In this paper, most depend ent variab les are normal- 
ized with respect to their boundary-layer edge values; p 
is scaled by p*u* 3 . Distances are scaled by the variable 
L* = (i/ e *x7 u*) 1/2 . The superscript * characterizes a di- 
mensional quantity, the subscript e indicates a value at the 
edge of a boundary layer, u is the streamwise velocity, and 
v is the kinematic viscosity. The viscosity, p* , is assumed 
to be given by the Sutherland formula. 


The techniques of compressible linear stability theory 
are quite well know; see, for example, [10]. The starting 
point is a laminar mean flow, q 0 . In most cases an approx- 
imate mean flow T is utilized, e.g., a quasi- parallel solution 
of the boundary-layer equations. One must then imagine 
that an appropriate forcing term has been added to Eq. 
(1). See [11] for a recent discussion. The total flow field, 
q(x,f), is written as 

q(x,<) = qo(*) + A(qi(r)e , '<°' +/3s '- u '‘ ) + c. c.) . (5) 

The streamwise and spanwise wave numbers are denoted 
by q and /?, respectively, and is the (temporal) fre- 
quency. The complex amplitude function, q x (z), deter- 
mines the structure in the wall-normal direction. The 
compressible Navier-Stokes equations are then linearized 
about the mean flow to first-order in the amplitude A. 
When combined with appropriate boundary conditions, an 
eigenvalue problem results. 

The spanvrise wavenumber is invariably taken to be real. 
In temporal theory, a is real and w is the (complex) eigen- 
value, with w the spatial growth rate. In the spatial con- 
cept, u is real and a is the eigenvalue, with -a, the spatial 
growth rate. In many cases the simpler temporal theory 
is applied, and the spatial growth rate approximated by 
— = Ui/Cg, where c g is the group velocity of the wave. 
The linear instability is referred to as the primary insta- 
bility. 

Some recent developments for compressible flow include 
incorporation of non-parallel effects through multiple-scale 
techniques [12], real gas effects [13], proper shock- wave 
boundary conditions [14], and clarification of some issues 
regarding propagation of three-dimensional waves [15]. 

The Af-factor method is applied by first computing the 
laminar mean flow past the body of interest and then 
applying LST to that flow. For a given frequency, w, 
the streamwise location at which a wave first becomes 
unstable, x 0 , is identified and then the spatial growth 
rate is integrated downstream to produce the A r -factor: 
A J(u) = J"* (-ai)dx'. (Note that if the amplitude of the 
instability at xo is denoted by Ao, then the amplitude at 
x is given by A/ Ac = thus, e N measures the am- 

plitude ratio.) This calculation is performed for a range of 
frequencies, and for each x, the maximum over u>, denoted 
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by just A\ is taken. 

The estimate of transition onset is based upon an em- 
pirical correlation between N and the location of transi- 
tion. The A r - factor method has had a surprising degree of 
success, even considering its limitations, such as inapplica- 
bility to flows in a high disturbance background for which 
the linear instability regime is “by-passed'’. The A r -factor 
method has matured to the point at which an analysis ca- 
pability is imminent for transition estimation across the 
speed range (including real gas effects) and for arbitrary 
steady three-dimensional mean flows. 

The computer requirements for solving a single LST 
eigenvalue problem are inconsequential. Even application 
of the A'-factor method to a three-dimensional mean flow 
requires only on the order of an hour of supercomputer 
time. However, a non-trivial related task is computa- 
tion of an accurate mean flow. LST requires accurate 
(and osciliation-free) mean flows and this is a far more 
stringent requirement than is customary in conventional 
steady-state CFD. This can take many tens of supercom- 
puter hours for a three-dimensional configuration. 

4 Secondary Instability Theory 


Secondary instability theory picks up where LST leaves 
off. In SIT one includes some weakly nonlinear (and three- 
dimensional) effects. One considers the linear stability 
with respect to secondary disturbances of a base flow com- 
prised of a laminar mean flow (assumed locally parallel) 
modulated by a small (but finite) amplitude primary dis- 
turbance. SIT is now well-established for incompressible 
flow. See Herbert [16] for a thorough review of the subject. 
In recent years it has been extended and applied to com- 
pressible boundary layers in [17], [18], [19], and [20]. Here 
we review some of the developments discussed in [20], 

The frequency, wavelength, and shape of the primary 
disturbance are obtained using LST. The primary wave 
is assumed to have no growth during the evolution of the 
secondary disturbance. In a Galilean frame, x + = x — c r <, 
moving with the phase velocity, c r , of the primary wave, 
the total flow variable, q, can be written in the Floquet 
form 

q=q ft + <3?{e"V'V‘ J =V' , ^ r+ f) q 2| ,(r)e^' r+ } , 

j = -oo 

( 6 ) 

where <\ h is given by Eq. (5). Hereafter, the subscripts 0, 
b, 1, and 2 refer, respectively, to the laminar mean flow, 
the modulated base flow, the primary disturbance, and the 
secondary disturbance. The primary wave amplitude, A , 



Figure 1: Variation of > with B when M e = 1.6. Funda- 
mental at A = 2.9%. Subharmonic at .4 = 1.5%. 

is defined so that it corresponds to the maximum value 
of the perturbation mass flux. The equations governing 
the secondary disturbance are obtained by linearizing Eq. 
(1) about the modulated base flow. The discretized dis- 
turbance equations are converted into a complex algebraic 
eigenvalue problem for cr or 7 and their associated eigen- 
functions. In temporal theory, 7 = 0 and cr ^ 0 is the 
complex eigenvalue to be determined. In the spatial con- 
cept, we write a — 7cy in Eq. (6) and solve for 7 as the 
eigenvalue [16]. The temporal and spatial growth rates are 
given by the real parts of cr and 7, respectively. The detun- 
ing parameter, h t defines the type of secondary instability. 
The subharmonic modes are given by h = 1, while the 
fundamental. modes correspond to h = 0. In practice the 
sum in Eq. (6) includes from 2 to 5 modes. 

As one example, consider the boundary layer over an 
insulated flat plate at Mach number M e = 1.6, Reynolds 
number Re = 1050, Prandtl number Pr = 0.72, and 
temperature T e = 216° Rankine. The primary wave 
is slightly damped, and is located near branch two of 
the neutral stability curve with a 1 = 0.1471 and F = 
10 6 x u\/Re = 82.6. For clarity, the secondary growth 
rate obtained from the temporal theory, which has been 
converted to a spatial growth rate by using the transfor- 
mation 7 r = cr r /c r , is termed the “transformed-spatial” 
growth rate. In many cases of interest both the dominant 
fundamental and subharmonic travel synchronously with 
the primary, i.e. , 7* = 0. 

In Fig. 1, the secondary growth rate is plotted 
-as a function of the normalized spanwise wavenumber 
B = 10 3 x (32 /Re. (For constant boundary-layer 

edge conditions, the parameter B signifies a wave of fixed 
spanwise wavelength as it propagates downstream.) The 
“transformed-spatial” growth rate agrees well with the 
(true) spatial growth rate. The most amplified subhar- 
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primary wave amplitude, A 


Figure 2: Variation of 7 r with A when A/ e — 1.6. 


monic disturbance consists of a pair of oblique waves trav- 
eling at equal and opposite angles at about 70° to the 
mean-flow direction. The most unstable fundamental dis- 
turbance is comprised of a stationary mode, which repre- 
sents a spanwise periodic distortion of the mean flow, and 
a pair of oblique waves propagating in opposite directions 
at about 57 ° to the mean-flow direction. 

Figure 2 depicts the typical catalytic effect of the pri- 
mary amplitude on the growth rates of the secondary dis- 
turbances. The spanwise wavenumbers, B , of the sub- 
harmonic and the fundamental are 0.19 and 0.22, respec- 
tively; these particular values of B correspond roughly to 
the most amplified secondary disturbances. The subhar- 
monic instability prevails over the fundamental instability, 
in an environment with a primary amplitude, .4, of less 
than about 2 . 8 %, while the converse is true for higher ’val- 
ues of -4. The transformed-spatial growth rates of the sub- 
harmonic modes are almost identital to the spatial rates 
for small -4 — the former increasingly underpredicts the 
latter as A increases from 1.5%. 

The preceding behavior is expected, since the use of 
temporal data to approximate spatial growth becomes less 
accurate when the growth rate is relatively high. Still, 
the transformation of temporal data to approximate the 
spatial growth rates of funda mental resonance modes is 
surprisinglv accurate even for a primary amplitude as high 
as 4%. 

The second example, given in Fig. 3, corresponds to a 
laminar boundary layer on an insulated 7-degree half-angle 
sharp cone at M t - 6 . 8 . The parameters are Re = 1939, 
T e = 128° Rankine, Pr = 0.70. The primary distur- 
bance is axisymmetric and is of the “second mode’’ type 
with o i = 0.2788 and F = 135. The two sets of cal- 
culations also serve to contrast the secondary instability 



Figure 3: Effect of the detuning parameter on 7 . M e = 6.8 
and B =0.135. 


arising from a first-mode primary (at A/ e — 1 . 6 ) with that 
arising from a second-mode primary (at A f e = 6 . 8 ). Over 
a range of Mach numbers up to at least 6 . 8 , subharmonic 
secondary instability ( h — 1 ) is found to prevail in a low 
primary-disturbance environment. This is especially true 
for high-speed flows. In particular, as illustrated in Fig. 
3 , no fundamental instability ( h = 0) for M e = 6.8 exists 
even for a primary amplitude of 2.625%. In fact, both 
the first- and second-mode primary waves associated with 
high-speed flows have been found ineffective in catalyzing 
unstable fundamental resonance modes. 

Recent developments in SIT include incorporation of 
non-parallel effects (but only for for the evolution of the 
primary wave) [17]. Ng & Erlebacher [19] have developed 
a fairly general capability which allows for oblique primary 
waves (important at low supersonic Mach number, where 
the most unstable primary is oblique). 

SIT has greater computational requirements than LST 
- the matrix eigenvalue problems which must be solved 
are larger than in LST. Nevertheless, a solution can be 
obtained in minutes on a supercomputer. 


5 Parabolized Stability Equations 


As a consequence of the rapid growth of the secondary 
wave, many harmonic waves, including the mean flow cor- 
rection mode, are excited to large amplitudes, and eventu- 
ally strongly nonlinear effects ensue; the flow then becomes 
transitional. Although SIT furnishes a guide to mecha- 
nisms that may be present near transition onset, it does 
not at present account for many non-parallel effects, and 
it incorporates only the lowest level non-linearity. 
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Figure 4: PSE computation of the effect of non-parallel 
mean flows for both 2-D and 3-D disturbances of a XI A == 
1.6 flow at F = 40. (Solid lines are from para! ■ - ^ T, 
dashed lines from PSE calculations, and symbols are irom 
multiple-scales analysis [12].) 


One efficient method which does incorporates both non- 
parallel and non-linear effects is the parabolized stabil- 
ity equations approach, first suggested by Herbert and 
Bertolotti [21], [22]. The PSE method facilitates the so- 
lution of the full partial differential equations for the dis- 
turbances by employing a partial parabolization along the 
dominant flow direction. In this approach, the disturbance 
is decomposed into a wave part and a shape function part. 
The elliptic terms are retained for the wave part, whereas 
the governing equations for the shape function are parab- 
olized in the stream wise direction. The parabolized equa- 
tions for the shape function are then solved numerically by 
a marching procedure. The wave properties are extracted 
from a local analysis. Nonlinear terms are formulated as 
forcing functions for the corresponding linear equations. 
Because the equation set contains non-parallel as well as 
nonlinear terms, the PSE method governs the spatial evo- 
lution of disturbances from the linear stage up to the tran- 
sitional stage. 

The PSE approach has been successfully applied to the 
stability of supersonic two-dimensional boundary layers by 
Bertolotti k Herbert [23] and Chang et ai [24]. In the 
linear regime, the method provides a means to include non- 
parallel effects due to the growth of the boundary layer, 
which is ignored in traditional LST. Furthermore, non- 



Figure 5: PSE amplitude evolution of the primary and 
subharmonic waves at AI C = 1.6, F = 50, and B = 0.053. 

linear effects such as wave/ wave interaction or secondary 
instability can be simulated by the non-linear PSE. 

To demonstrate the applicability of the PSE approach, 
we present some of the results given in [24]. Linear 
PSE calculations were performed for a Mach 1.6 flat-plate 
boundary layer previously studied by EI-Hady [12]. The 
frequency of the disturbances is F = 40. Calculations 
w r ere performed for both 2-D and 3-D linear disturbances; 
the wave angle for the oblique, 3-D waves was about 50°. 
The growth rate of the mass flow fluctuations from the 
PSE calculations together with the multiple-scales results 
are plotted along with the growth rates obtained by quasi- 
parallel LST in Fig. 4. PSE results agree quite well with 
those obtained from the multiple-scales approach. The re- 
sults also indicate that for the first mode disturbance at 
Mach 1.6, flow non-parallelism has more effect on three- 
dimensional disturbances than on two-dimensional ones. 
The non-parallel effect on oblique waves is qualitatively 
very similar to that in incompressible flows, as found by 
Bertolotti [21]. 

Compressible non-linear PSE computations for sec- 
ondary instability mechanisms and the subsequent start 
of laminar breakdown have also been demonstrated . The 
flow is again a Mach 1.6 flat-plate boundary layer with a 
primary disturbance frequency of F = 50. The free-stream 
temperature is 540 ° Rankine and the Prandtl number is 
Pr = 0 . 71 . Figure 5 shows the evolution of primary and 
subharmonic disturbances for various initial amplitudes of 
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Figure 6: PSE amplitude evolution of various Fourier 

modes for a M e = 1.6 fundamental breakdown 


the primary waves (the subharmonic amplitudes are the 
same for all three cases). The spamvise wave number of 
the subharmonic mode is fixed at B = 0.053, which corre- 
sponds to an oblique wave angle of 45°. As can be seen, a 
1.1% initial amplitude for the primary mode is enough to 
trigger the secondary growth. 

Non-linear PSE calculations were also performed for the 
same Mach 1.6 case but for a fundamental-type secondary 
resonance. The initial amplitude of the primary wave is 
again 3% and that of the secondary is taken to be 0.005%- 
The spanwise wave number is B — .152 (oblique wave an- 
gle of 60° for the secondary wave) and the primary wave 
frequency is again F — 50. The non-linear evolution of 
the maximum rms amplitude of u ; (a prime is used to de- 
note the fluctuating part of a variable) is shown in Figure 
6. Clearly, the presence of the primary 2-D and secondary 
3-D disturbances results in wave resonance and strong sec- 
ondary growth of the oblique wave. When the secondary 
disturbance is amplified to about the same amplitude of 
the primary wave, all harmonics are excited and the flow 
becomes transitional. This is confirmed by plotting the 
average wall shear in Figure 7. The computed wall shear 
is only slightly above the laminar value for most of the 
computational domain. (The PSE wall shear lies above 
the laminar value right from the beginning because of the 
relatively high initial amplitude of the 2-D primary dis- 
turbance.) Eventually the wall shear departs sharply from 
the laminar value, indicative of transition onset. Thus, the 
PSE computation captures the skin friction rise, which is 
one of the criteria for transition onset) accurate predic- 
tion of its location is a prime goal of transition prediction 
methods. 
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Figure 7: PSE wall shear of laminar and perturbed flows 
for a M e = 1.6 fundamental breakdown. 

Linear PSE computations (equivalent to non-parallel 
LST) are quite cheap. Nonlinear PSE can compute up 
to the skin friction rise in no more than an hour of su- 
percomputer time. However, current numerical techniques 
for nonlinear PSE have computational requirements which 
scale quadratically with both the number of spanwise 
Fourier components and the number of temporal frequen- 
cies retained in the approximation. This means that PSE 
computations for the later stages of transition and for ran- 
dom inflow and/or freestream conditions are exceedingly 
expensive. 

The SIT and PSE tools that have been described up 
to this point are oriented towards forced transitional. e., 
transition characterized by sharply defined frequencies as 
might occur from specific forcing. The technologically in- 
teresting problem is that of natural transition, for which 
a broad range of frequencies is present. To capture the 
nonlinear interaction between a wide range of frequencies, 
DNS and LES are currently the most appropriate tools. 


6 Direct Numerical Simulation 

Direct numerical simulation solves the time-dependent, 
three-dimensional, nonlinear, Navier-Stokes equations 
subject to prescribed initial and boundary conditions with- 
out recourse to empirical models. A thorough review of 
this approach has been given by Kleiser & Zang [25]. 
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Figure 8: Spatial SIT vs. DNS for a M e = 1.6 fundamental 
instability. 

When attacked as the total bound ary- layer transition 
problem - from receptivity through fully-developed tur- 
bulence - the non-parallel, spatial formulation is certainly 
more appropriate than the temporal approximation. How- 
ever, DNS for the complete transition process is an exceed- 
ingly expensive tool even for the low Reynolds numbers 
to which it is of necessity restricted. Gilbert k Kleiser 
[26] performed the first well-resolved simulation of the 
complete transition to turbulence and this took several 
hundred supercomputer hours for a temporal computa- 
tion of forced incompressible transition. With somewhat 
relaxed resolution requirements Rai k Moin [27] have re- 
cently computed bypass transition for low-speed flow past 
a flat plate. This required many hundreds of supercom- 
puter hours and it remains to be seen what the require- 
ments are for a well-resolved computation for this problem. 
For the foreseeable future, both temporal and spatial DNS 
have a role, but this tool ought to be applied selectively. 

One role for DNS is the corroboration of simpler tools, 
such as SIT and PSE. For compressible flows it has been 
used to verify temporal SIT [19], [9], spatial SIT [20], and 
some aspects of nonlinear PSE [24]. Given the algebraic 
complexity of SIT and PSE, particularly for oblique pri- 
maries, this role is a needed one to establish confidence in 
them. (It also furnishes a stringent calibration of DNS ) 

One comparison between spatial (but quasi-parallel) 
SIT and DNS by Ng k Zang [20] was performed for a 
fundamental type instability at M c = 1.6, Re = 613, 
Pr = 0.70, and T t = 520° Rankine. The primary wave is a 
2-D first mode with frequency F = 73. The subharmonic 
secondary wave has spanwise wavenumber B — 0.1465. 
The amplitude of the primary is chosen to be 6%, while 
that of the secondary is 0.1%. Although the spatial DNS 
code is designed for non-parallel flow, for comparison with 



Figure 9: Evolution of selected harmonics from DNS of 
M e =4.5 cylinder subharmonic transition. 


the quasi-parallel SIT theory the mean flow was con- 
strained to be parallel by using a forcing function. The ini- 
tial conditions consisted solely of the laminar mean flow. 
At the inflow boundary, x m = Xq = uJL* /*£, the flow 
consisted of the mean flow plus the appropriate contribu- 
tions from the linear and secondary eigenfunctions. The 
physical domain consisted of 8 wavelengths of the pri- 
mary wave, with ‘'buffer domain” modifications [28] to the 
Navier-Stokes equations in the last 2 wavelengths used to 
ameliorate potential difficulties with the outflow boundary 
conditions. (See [20] for details.) Figure 8 shows the span- 
wise velocity component, v, at z = 0.26*, where 6* is the 
displacement thickness, and fay = ?r/2 after 10 periods of 
forcing; this component is due solely to the secondary in- 
stability. The agreement between the DNS and SIT results 
is excellent, except, of course, in the buffer domain. 

This computation utilized 12 points per stream wise 
wavelength. For transition in high-speed flows the growth 
rates of both primary and secondary disturbances are 
lower than for incompressible flow. As a consequence on 
the order of 10 2 wavelengths would be needed to follow the 
prim ary /secondary stages from about the 1% level to the 
skin friction rise. This is a prohibitive expense and argues 
strongly for the use of simpler methods such as SIT and 
especially PSE for routine application to the early stages 
of transition. 

The unique niche for DNS is the highly nonlinear, lami- 
nar breakdown stage and the subsequent transition to tur- 
bulence; for this the non-parallel effects appear to be less 
significant than they are for the rather lengthy primary 
and secondary instability stages. These early stages are 
nowadays treated far more efficiently by SIT and PSE ap- 
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Figure IQ: Skin friction and shape factor from DNS of 
M e = 4.5 cylinder subharmonic transition. 


proaches than by DNS. Indeed, the simpler theories can 
well be used to set the stage for the DNS. LST selects the 
dominant primary instability (and determines the relevant 
streamwise scales); SIT and/or PSE select the dominant 
spanwise scales and can be used to “jump start 5 the DNS 
at fairly high disturbance levels. 


This is the approach that was taken by Pruett Zang 
[29] in their temporal DNS of transition in Mach 4.5 flow 
past a cylinder. The primary disturbance was a second 
mode and the secondary disturbance was of subharmonic 
type. Due to the periodicity assumptions in the stream- 
wise (j:) and spanwise (y) directions, the dependent vari- 
ables have Fourier series representations in these direc- 
tions. A useful measure of the strength of a given Fourier 
harmonic is 


aOO 

E kl ,k v (t) = J />o(=)|u*„*,M)| : dz, (7) 


where Uk x ,k y is the Fourier coefficient of the velocity cor- 
responding to wavenumbers k x and k y (with respect to ai 
and /?2 1 respectively; 

4 = 2 - ho ( 8 ) 

accounts for some of the symmetries in the problem. The 
quantity Ek X} k y is approximately the kinetic energy of the 
(&*, k y ) mode. 

Figure 9, taken from [29], summarizes the time evolu- 
tion of the principal modes for the Mach 4.5 cylinder sim- 



Figure 11: Cartoon of the precursor transition effect (top) 
[30] and its manifestation in the spatially-reconstructed 
Reynolds stress from DNS of M t - 6.8 cone subharmonic 
transition (bottom). 


ulation. (Time is reckoned in units of the period of the 
primary wave.) The DNS proceeded through the stages 
of primary /secondary instability then underwent weakly 
nonlinear and strongly nonlinear stages, and finally com- 
menced laminar breakdown. The stages cited above last 
from 0-15, 15-35, 35-45 and 45-60, respectively. The sym- 
bols on the figure are the predictions from LST and SIT 
for the growth of the primary and secondary disturbances. 
They are in good agreement with the DNS. One interest- 
ing feature of this simulation is the prominence that the 
(0, 2) mode assumes in the latter stages of transition. This 
mode is not present in the initial conditions (nor in SIT) 
and is generated by nonlinear effects. Additional DNS 
have revealed that this mode plays an essential role in the 
final laminar breakdown. 

Figure 10 shows the evolution of the skin friction, C/, 
and the shape factor, H , for the Mach 4.5 cylinder tran- 
sition. This simulation was stopped at about 60 peri- 
ods because of strong gradients that even its million grid 
points could not resolve. Indeed, the judgment of Pruett 
& Zang [29] was that the resolution became questionable 
after 55 periods. This resolution problem intensifies as 
Mach number increases, and in a particular computation 
may eventually manifest itself in negative values of some 
of the thermodynamic quantities. This particular diffi- 
culty does not arise for low-speed DNS. For compressible 
flow the dilemma is that highly- accurate central-difference 
schemes do not have sufficient artificial viscosity to resolve 
strong gradients at high Reynolds number, whereas con- 
ventional upwind CFD schemes are so dissipative that they 
corrupt the delicate physics of transition. An encourag- 
ing recent development is the work of Atkins [30], which 
demonstrated good results for a compressible free shear 
layer transition using a fifth-order ENO scheme. 

The shape factor plot suggests that the simulation 
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Figure 12: Experimental (top) [32] and DNS (bottom) 
[25] visualizations of the rope-like structures in supersonic 
transition. 


has proceeded almost completely through the transition 
process. Figure 11 presents the spatially-reconstructed 
streamwise velocity fluctuations, as represented by the 
Reynolds-stress component Tu = — ~pu / u r , from a Mach 
6.8 cone simulation [29]. These exhibit the so-called “pre- 
cursor transition effect”, sketched in the cartoon in the top 
part of the figure, whereby the transition originates near 
the boundary-layer edge and gradually propagates toward 
the wall. This same effect is present in the DNS. Flow- 
field visualizations presented in Fig. 12 demonstrate the 
presence in the Mach 4.5 cylinder DNS of the “rope-like 
structures” that have been observed in numerous experi- 
ments ([31], [32], [33]). One of the more significant conclu- 
sions of [29] was that the rope-like structures are actually 
manifestations of SIT and not LST, as had long been sus- 
pected. 

A last sample of results from this DNS is presented in 
Fig. 13. It illustrates the evolution of the turbulent kinetic 
energy, k , and the turbulent dissipation, e, through the 
transition region. This kind of information, supplemented 
by detailed information on the key terms in the evolution 
equations for these quantities, has the potential to lead to 
substantial improvements in two-equation models for tran- 
sitional flow. However, DNS (particularly spatial DNS) is 
so computationally intensive that a less drastic, but still 
effective, tool for exploring the later stages of transition 
would be quite desirable. The following section describes 



0 12 3 

z'/6' 


Figure 13: Evolution of turbulent kinetic energy and dis- 
sipation for M e = 4.5 cylinder subharmonic transition. 

a promising candidate. 


7 Large-Eddy Simulation 

In large-eddy simulation the small scales of the flow are 
modeled in terms of the large-scale flow. The Reynolds- 
averaged Navier-Stokes equations, in contrast, model the 
higher-order moments in terms of the lower-order mo- 
ments. In LES the flow variables are decomposed into 
a large-scale (resolvable) component and a small-scale 
(subgrid-scale) component. LES was originally developed 
for turbulent flow (see [34] for a survey of the state-of- 
the-art in LES), and, at least with the more established 
subgrid-scale (SGS) models, some refinements have proven 
necessary to handle transitional flow properly. In the tran- 
sitional case one not only wants to have a model which 
works well for the final turbulent state, but also one which 
captures the primary, secondary and nonlinear stages cor- 
rectly (without, for example, exerting excessive damping 
of the instability waves), predicts well the location of tran- 
sition onset, and makes good predictions from transition 
onset through the transitional zone to the fully turbulent 
state. Piomelli and co-workers ([35], [36], [37], [38]) have 
led the effort to utilize DNS of transition to calibrate and 
refine SGS models for this process. This work has to date 
been confined to incompressible flow. In this section we 
describe some recent developments for compressible tran- 
sition modeling via LES. 
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The large-scale field is defined by the filtering operation: 


JF(x) = J C(x,x’)T (x')fx', 


where the integral is extended over the entire spatial do- 
main and G = GiG 2 G 3 , where G,(x,) is the filter func- 
tion in the ith direction. For the velocity and tempera- 
ture, Favre filtering is utilized: T = f + T ' , where T' 

is the SGS part of and the Favre filter is defined by 

f - JiF/~p. The sharp Fourier cutoff filter is chosen for 
this work because of previous experience with this filter in 
LES of incompressible transition ([35]). 

The dimensionless governing equations for compressible 

LES are 

8rp djfiuk) _ 

dt dxk 

d(puh) djpukui) <9p dau dni 

dt dxi ~ dx k dxi dxi 

dp _ dp_ _duk 1 dh 

dt + Uk dx k + 7P dx k M; PrRedx k 

+ lzlf+ i ?2i + (T - 1)3. (, - 1)8, ^2- 

Re M} c)x k ox k ox k 

yMeP = pT . (10) 


S k i = (duk/dxi +dui/dxk)J 2 is the Fav re-filtered rate- 
of-strain tensor and II § = S r nnS mn is its second invariant. 
Cr is the compressible Smagorinsky constant, Prj is the 
turbulent Prandtl number, and A = (A.r AyAr) 1 / 3 . The 
function T is an intermittency-like term that turns itself 
on slowly in the transitional zone [35]. Eq. (13) represents 
a Van Driest wall damping and z + indicates a wall-normal 
distance made dimensionless by the wall shear velocity and 
kinematic viscosity. 

The second SGS model considered here is the structure 
function model [43], which is based on a physical space 
implementation of the concept of spectral eddy viscosity. 
Some results for this model have been reported for a spa- 
tial LES of a Mach 5 boundary layer [41]. The structure 
function model is of the following form for the SGS shear 
stress and heat flux: 

n , = CrVT p (x, Ax, Ay, t) x f§ki - ^ Srnm Ski 

(14) 


Qk = — 


CrVTpAx /= 


F 2 (x. Ax, Ay, Gtt 


The SGS stress tensor r k i and the j>GS heat flu x Q k are 
defined by r k \ = — p(u*ur- UfcWf + u k uj + uju* + uj.u{) and 

Qk — —J>{ UkT — u k T + u f k T T- u k T f + tt / fc T / ). 

There have been a number of SGS models proposed for 
compressible LES ([39], [40], [41], [42]). Two of these mod- 
els have been applied to the Mach 4.5 transition problem 
discussed in the preceeding section. 

The first SGS model considered is the SEZHu model 
derived by Speziale et al. [40] for compressible isotropic 
turbulence. This model has been chosen because there 
are now available some extensive a posteriori comparisons 
of its performance on compressible, isotropic turbulence 
([44], [45]). Following the work of Piomelli et al. [46], only 
the Smagorinsky portion of the SEZHu model is used with 
the Fourier cutoff filter. (This filter is applied in x and y; 
no filtering is applied in the inhomogeneous z direction.) 
Hence, the SGS stress model is of the form 

r kl = 2C R VTpA 2 Il\ /2 (§ki - , ( 11 ) 


and the SGS heat flux is given by 



Qt =- c f T *A’„yf, 

^ Pr T s dxk 

(12) 

where 




V = (l -e" z+/25 ) , 

(13) 


where 

F 2 = 5 ([]u(jj -F Ax,y,z,<)-u(x,y,z,0ll 2 

+ ||u(x,y, z,t) - \i(x - &x,y, z,t)\\ 2 

+ ^(||u(x,y + Ay, z,t) - u (x,y, z , t )|| 2 

Ay 

+ ||u(x,y, z,i) — u(x,y- Ay, z,<)|| 2 )) . (16) 

Although neither Van Driest wall damping nor the inter- 
mittency term were in the structure function model as 
used in [41], they were added here as they seem to furnish 
better results. 

The structure-function SGS model was tested a poste- 
riori both with and without the “intermittency” function^ 
and a comparison is given in Pig. 14 for the primary and 
secondary components. (In the latter case T = 1.) The 
coefficients for these runs were Cr = 0.06, which is the 
value recommended in [41], and Prr = 0.70. The LES 
used IQ 4 grid points, two orders of magnitude fewer than 
the DNS. The original structure function model is clearly 
far too dissipative in the early stages, whereas the modi- 
fied model agrees very well with the high- resolution DNS. 
In this respect these results are similar to those reported 
in [35] for the original and modified Smagorinsky model 
when applied to incompressible transition. The integral 
properties are in quite good agreement with those of the 
high-resolution DNS up to T = 55; they are summarized 
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Figure 14: Comparison of unmodified and modified struc- 
ture function for LES of M e = 4.5 cylinder transition. 


in Figure 15. Note that the computation reaches the fully 
turbulent regime. 

For the SEZHu model the nominal coefficient is Cr = 
0.011 [44]. However, on the same grid as was used for 
the previous LES computations this had to be increased 
to Cr = 0.50 in order to compute all the way through 
to turbulence. The reason for this is not entirely clear at 
present. It might be that the larger value of the Smagorin- 
sky constant serves to provide artificial viscosity needed to 
stabilize the computation, or it might be due to the much 
smaller length scales involved in second mode transition. 
The later stages of the DNS suggest that for the turbulent 
state the computational box, in wall units, was 240 in x 
and 150 in y ) as opposed to typical values from incom- 
pressible flow of 2000 and 1000, respectively. 

Krai & Zang [47] have performed some LES of a Mach 
4.5 turbulent boundary layer with computational domains 
closer to the standard incompressible sizes. Here they 
found that reasonable results could be obtained with the 
constants Cr closer to the accepted incompressible values. 
It appears that application of the dynamic eddy viscosity 
concept [37] to this problem would be quite fruitful. 

The potential of LES for transition is that it permits 
computations through the transitional zone at an order of 
magnitude or more lesser expense than for DNS. A discus- 
sion of what sort of information can be reliably provided 
by LES and DNS is given in [48]. For incompressible flow, 
demonstrations are needed for spatial transition; for com- 




Figure 15: Evolution of skin friction and shape factor for 
LES of M e = 4.5 cylinder transition. 

pressible flow, the role of the SGS viscosity needs clarifi- 
cation: To what extent is it furnishing artificial viscosity 
rather than serving its intended purpose of modeling the 
physics? 


8 Prospects 


The past several years have witnessed many exciting de- 
velopments in transition research, particularly for com- 
pressible flow, that make an absolute amplitude criterion 
for prediction of transition onset a tantalizing prospect. 
Many of the components of such a methodology have been 
discussed in this paper. LST technology is virtually com- 
plete for real configurations. SIT is likely to be absorbed 
within PSE. The scope of PSE needs to be vastly increased 
and it would greatly benefit from a firmer mathemati- 
cal foundation. DNS will no doubt undergo algorithmic 
improvements and take advantage of increased computa- 
tional power. LES will likely evolve through several gen- 
erations of SGS model improvements. 

The philosophy, not only for transition prediction, but 
also for basic research into transition physics, ought to be 
to apply to each stage of transition or to each physical 
problem the most economical and revealing method in the 
transition prediction toolkit. 
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